#include <iostream>
#include <vector>
#include <cmath>
#include "matplotlibcpp.h"

namespace plt = matplotlibcpp;

// 定义截断幂函数
double truncated_power(double x, int n, double t) {
    if (x > t) return 0.0;
    return std::pow(t-x, n);
}

// 计算除差
double divided_difference_1(double x, int n, double t1, double t2) {
    double diff;
    diff = (truncated_power(x, n, t2) - truncated_power(x, n, t1))/(t2 - t1);
    return diff;
}

// 计算除差
double divided_difference_2(double x, int n, double t1, double t2, double t3) {
    double diff;
    diff = (divided_difference_1(x, n, t2, t3) - divided_difference_1(x, n, t1, t2))/(t3 - t1);
    return diff;
}

// 计算除差
double divided_difference_3(double x, int n, double t1, double t2, double t3, double t4) {
    double diff;
    diff = (divided_difference_2(x, n, t2, t3, t4) - divided_difference_2(x, n, t1, t2, t3))/(t4 - t1);
    return diff;
}

void plot11(){
    // 定义节点
    std::vector<double> x;
    std::vector<double> y1, y2, y3;

    // 计算n=1和n=2的截断幂函数值
    int N = 100;
    double h = 3.0/N;
    double ti_1 = 1, ti = 2, ti1 = 3;
    for (int i = 0; i<=N; i++) {
        double xi = i*h;
        x.push_back(xi);
        y1.push_back(truncated_power(xi, 1, ti_1));
        y2.push_back(truncated_power(xi, 1, ti));
        y3.push_back(truncated_power(xi, 1, ti1));
    }

    // 绘制图形
    
    plt::figure_size(600, 200);
    plt::plot(x, y1, "r-");
    plt::xlim(0, 3);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_1_11.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y2, "r-");
    plt::xlim(0, 3);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_1_12.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y3, "r-");
    plt::xlim(0, 3);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_1_13.png");
    plt::show();
}

void plot12(){
    // 定义节点
    std::vector<double> x;
    std::vector<double> y1, y2;

    // 计算n=1和n=2的截断幂函数值
    int N = 100;
    double h = 3.0/N;
    double ti_1 = 1, ti = 2, ti1 = 3;
    for (int i = 0; i<=N; i++) {
        double xi = i*h;
        x.push_back(xi);
        y1.push_back(divided_difference_1(xi, 1, ti_1, ti));
        y2.push_back(divided_difference_1(xi, 1, ti, ti1));
    }

    // 绘制图形
    
    plt::figure_size(600, 200);
    plt::plot(x, y1, "r-");
    plt::xlim(0, 3);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_1_21.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y2, "r-");
    plt::xlim(0, 3);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_1_22.png");
    plt::show();

}

void plot13(){
    // 定义节点
    std::vector<double> x;
    std::vector<double> y1;

    // 计算n=1和n=2的截断幂函数值
    int N = 100;
    double h = 3.0/N;
    double ti_1 = 1, ti = 2, ti1 = 3;
    for (int i = 0; i<=N; i++) {
        double xi = i*h;
        x.push_back(xi);
        y1.push_back(divided_difference_2(xi, 1, ti_1, ti, ti1));
    }

    // 绘制图形
    plt::figure_size(600, 200);
    plt::plot(x, y1, "r-");
    plt::xlim(0, 3);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_1_31.png");
    plt::show();
}

void plot21(){
    // 定义节点
    std::vector<double> x;
    std::vector<double> y1, y2, y3, y4;

    // 计算n=1和n=2的截断幂函数值
    int N = 100;
    double h = 4.0/N;
    double ti_1 = 1, ti = 2, ti1 = 3, ti2 = 4;
    for (int i = 0; i<=N; i++) {
        double xi = i*h;
        x.push_back(xi);
        y1.push_back(truncated_power(xi, 2, ti_1));
        y2.push_back(truncated_power(xi, 2, ti));
        y3.push_back(truncated_power(xi, 2, ti1));
        y4.push_back(truncated_power(xi, 2, ti2));
    }

    // 绘制图形
    
    plt::figure_size(600, 200);
    plt::plot(x, y1, "r-");
    plt::xlim(0, 4);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_11.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y2, "r-");
    plt::xlim(0, 4);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_12.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y3, "r-");
    plt::xlim(0, 4);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_13.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y4, "r-");
    plt::xlim(0, 4);
    plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_14.png");
    plt::show();
}

void plot22(){
    // 定义节点
    std::vector<double> x;
    std::vector<double> y1, y2, y3;

    // 计算n=1和n=2的截断幂函数值
    int N = 100;
    double h = 4.0/N;
    double ti_1 = 1, ti = 2, ti1 = 3, ti2 = 4;
    for (int i = 0; i<=N; i++) {
        double xi = i*h;
        x.push_back(xi);
        y1.push_back(divided_difference_1(xi, 2, ti_1, ti));
        y2.push_back(divided_difference_1(xi, 2, ti, ti1));
        y3.push_back(divided_difference_1(xi, 2, ti1, ti2));
    }

    // 绘制图形
    
    plt::figure_size(600, 200);
    plt::plot(x, y1, "r-");
    plt::xlim(0, 4);
    // plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_21.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y2, "r-");
    plt::xlim(0, 4);
    // plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_22.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y3, "r-");
    plt::xlim(0, 4);
    // plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_23.png");
    plt::show();
}

void plot23(){
    // 定义节点
    std::vector<double> x;
    std::vector<double> y1, y2;

    // 计算n=1和n=2的截断幂函数值
    int N = 100;
    double h = 4.0/N;
    double ti_1 = 1, ti = 2, ti1 = 3, ti2 = 4;
    for (int i = 0; i<=N; i++) {
        double xi = i*h;
        x.push_back(xi);
        y1.push_back(divided_difference_2(xi, 2, ti_1, ti, ti1));
        y2.push_back(divided_difference_2(xi, 2, ti, ti1, ti2));
    }

    // 绘制图形
    
    plt::figure_size(600, 200);
    plt::plot(x, y1, "r-");
    plt::xlim(0, 4);
    // plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_31.png");
    plt::show();

    plt::figure_size(600, 200);
    plt::plot(x, y2, "r-");
    plt::xlim(0, 4);
    // plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_32.png");
    plt::show();

}

void plot24(){
    // 定义节点
    std::vector<double> x;
    std::vector<double> y1;

    // 计算n=1和n=2的截断幂函数值
    int N = 100;
    double h = 4.0/N;
    double ti_1 = 1, ti = 2, ti1 = 3, ti2 = 4;
    for (int i = 0; i<=N; i++) {
        double xi = i*h;
        x.push_back(xi);
        y1.push_back(divided_difference_3(xi, 2, ti_1, ti, ti1, ti2));
    }

    // 绘制图形
    
    plt::figure_size(600, 200);
    plt::plot(x, y1, "r-");
    plt::xlim(0, 4);
    // plt::ylim(0, 1);
    plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/PF_2_41.png");
    plt::show();

}

int main() {
    plot11();
    plot12();
    plot13();
    plot21();
    plot22();
    plot23();
    plot24();
    // plt::figure_size(600, 200); // 设置图形的尺寸
    // plt::save("F:/Codefield/numerical_analysis/Programming/Chapter3/figure/blank.png");
    // plt::show();
    return 0;
}